res=GET('https://developer.nps.gov/api/v1/parks?limit=500&api_key=B9nDpbkbrb3kSOjz6kXSxMJ3d6MSpUvt1QqYdeyn')

data = res %>% content("text") %>% jsonlite::fromJSON() %>% as_tibble()

NPS_data=data %>% unnest(data) %>% select(fullName,latitude,longitude,activities,states, parkCode) %>%  janitor::clean_names() %>% 
  mutate(
    latitude = as.numeric(latitude), 
    longitude = as.numeric(longitude)
  )
visitation_data <- 
  read_csv("data/Query Builder for Public Use Statistics (1979 - Last Calendar Year).csv") %>% 
  janitor::clean_names() %>% 
  mutate(unit_code = tolower(unit_code)) %>% 
  rename(park_code = unit_code, 
         full_name = park_name) %>% 
  select(full_name, park_code, park_type, region, state, year, month, recreation_visits, tent_campers, rv_campers, tent_campers, backcountry)
## Rows: 101419 Columns: 35
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (10): ParkName, UnitCode, ParkType, Region, State, ParkNameTotal, UnitCo...
## dbl  (3): Year, Month, YearTotal
## num (22): RecreationVisits, NonRecreationVisits, RecreationHours, NonRecreat...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
combined_data <- right_join(NPS_data, visitation_data, by = c("park_code"))

faceted plot with rv, tent, and backcountry visits

visitation_data %>% 
  group_by(park_type, month) %>% 
  summarize(total_visitation = sum(recreation_visits)) %>% 
plot_ly(x = ~as.factor(month), y = ~total_visitation, type = "scatter", mode = "lines", color = ~park_type)
## `summarise()` has grouped output by 'park_type'. You can override using the
## `.groups` argument.
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors
visit_summary <- 
  visitation_data %>% 
  group_by(park_type, month) %>% 
  summarize(total_tent = sum(tent_campers), 
            total_backcountry = sum(backcountry), 
            total_rv = sum(rv_campers))
## `summarise()` has grouped output by 'park_type'. You can override using the
## `.groups` argument.
tent_campers_plot <- plot_ly(visit_summary, 
                             x = ~as.factor(month), 
                             y = ~total_tent, 
                             type = "scatter", 
                             mode = "lines", 
                             color = ~park_type,
                             legendgroup = ~park_type) %>%
  layout(title = "Tent Campers")

backcountry_plot <- plot_ly(visit_summary, 
                            x = ~as.factor(month), 
                            y = ~total_backcountry, 
                            type = "scatter",
                            mode = "lines",
                            color = ~park_type,
                           legendgroup = ~park_type,
                           showlegend = FALSE) %>%
  layout(title = "Backcountry")

rv_campers_plot <- plot_ly(visit_summary, 
                           x = ~as.factor(month), 
                           y = ~total_rv, 
                           type = "scatter", 
                           mode = "lines",
                           color = ~park_type,
                          legendgroup = ~park_type,
                          showlegend = FALSE) %>%
  layout(title = "RV Campers")

subplot(
  tent_campers_plot, 
  backcountry_plot, 
  rv_campers_plot, 
  nrows = 2,  
  shareX = TRUE, 
  shareY = TRUE   
)
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors

visits by season boxplot

visitation_data %>% 
  mutate(season = case_when(
    month %in% c(12,1,2) ~ "Winter", 
    month %in% c(3,4,5) ~ "Spring", 
    month %in% c(6,7,8) ~ "Summer", 
    TRUE ~ "Fall"
  )) %>% 
  plot_ly(x = ~season, y = ~recreation_visits, type = "box")

total visitation by year

visitation_data %>% 
  group_by(park_type, year) %>% 
  summarize(total_visitation = sum(recreation_visits)) %>% 
  ggplot(aes(x = as.factor(year), y = total_visitation, color = park_type, group = park_type)) + 
  geom_line(linewidth = 0.7) +  theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1))
## `summarise()` has grouped output by 'park_type'. You can override using the
## `.groups` argument.

total visitation by season

visitation_data %>% 
  mutate(season = case_when(
    month %in% c(12,1,2) ~ "Winter", 
    month %in% c(3,4,5) ~ "Spring", 
    month %in% c(6,7,8) ~ "Summer", 
    TRUE ~ "Fall"
  )) %>% 
  group_by(season) %>% 
  summarize(total_tent = sum(tent_campers), 
            backcountry_visits = sum(backcountry), 
            total_rv = sum(rv_campers)) %>% 
  pivot_longer(
    total_tent:total_rv, 
    values_to = "total_visit", 
    names_to = "type_visit",
    names_prefix = "total_"
  ) %>% 
  plot_ly(x = ~season, y = ~total_visit, color = ~type_visit, type = "bar")

mean visitation by park type and season

visitation_data %>% 
  mutate(season = case_when(
    month %in% c(12,1,2) ~ "Winter", 
    month %in% c(3,4,5) ~ "Spring", 
    month %in% c(6,7,8) ~ "Summer", 
    TRUE ~ "Fall"
  )) %>% 
  group_by(park_type, season) %>% 
  mutate(
    mean_total = mean(recreation_visits, na.rm = TRUE)
  ) %>% 
    plot_ly(x = ~season, y = ~mean_total, color = ~park_type)
## No trace type specified:
##   Based on info supplied, a 'bar' trace seems appropriate.
##   Read more about this trace type -> https://plotly.com/r/reference/#bar
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors
## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors

mean visitation season and visit type

visitation_data %>% 
  mutate(season = case_when(
    month %in% c(12,1,2) ~ "Winter", 
    month %in% c(3,4,5) ~ "Spring", 
    month %in% c(6,7,8) ~ "Summer", 
    TRUE ~ "Fall"
  )) %>% 
  group_by(season) %>% 
  summarize(mean_tent = mean(tent_campers), 
            mean_backcountry = mean(backcountry), 
            mean_rv = mean(rv_campers)) %>% 
  pivot_longer(
    mean_tent:mean_rv, 
    names_to = "type_visit", 
    values_to = "mean", 
    names_prefix = "mean_"
  ) %>% 
  plot_ly(x = ~season, y = ~mean, color = ~type_visit)
## No trace type specified:
##   Based on info supplied, a 'bar' trace seems appropriate.
##   Read more about this trace type -> https://plotly.com/r/reference/#bar

map of all national parks in US

leaflet(data = NPS_data) %>% 
  addTiles() %>% 
  addMarkers(~longitude, ~latitude) %>%
  setView(lng = -98.583, lat = 39.828, zoom = 4)

Mean recreation visit by park type (separate plots for each park type)

park_types <- unique(visitation_data$park_type)
plots <- list()
for (i in seq_along(park_types)) {
plots[[i]] <- visitation_data %>% 
    filter(park_type == park_types[i]) %>% 
    mutate(season = case_when(
      month %in% c(12, 1, 2) ~ "Winter", 
      month %in% c(3, 4, 5) ~ "Spring", 
      month %in% c(6, 7, 8) ~ "Summer", 
      TRUE ~ "Fall"
    )) %>% 
    group_by(season) %>%
    summarize(avg_visits = mean(recreation_visits, na.rm = TRUE)) %>%
    ggplot(aes(x = season, y = avg_visits, fill = season)) +
    geom_col() +
    ggtitle(paste("Average Visits by Season for", park_types[i])) 
} 
plots
## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

## 
## [[7]]

## 
## [[8]]

## 
## [[9]]

## 
## [[10]]

## 
## [[11]]

## 
## [[12]]

## 
## [[13]]

## 
## [[14]]

## 
## [[15]]

## 
## [[16]]

## 
## [[17]]

## 
## [[18]]

National Preserves has pretty similar avg visits across seasons

Average visitation by type and park_type (separate graph for each park_type)

park_types <- unique(visitation_data$park_type)
plots <- list()

for (i in seq_along(park_types)) {
  visit_summary <- visitation_data %>% 
    filter(park_type == park_types[i]) %>% 
    mutate(season = case_when(
      month %in% c(12, 1, 2) ~ "Winter", 
      month %in% c(3, 4, 5) ~ "Spring", 
      month %in% c(6, 7, 8) ~ "Summer", 
      TRUE ~ "Fall"
    )) %>% 
    group_by(season) %>% 
    summarize(
      mean_tent = mean(tent_campers, na.rm = TRUE), 
      mean_backcountry = mean(backcountry, na.rm = TRUE), 
      mean_rv = mean(rv_campers, na.rm = TRUE)
    ) %>% 
    pivot_longer(
      cols = starts_with("mean_"), 
      names_to = "type_visit", 
      values_to = "mean", 
      names_prefix = "mean_"
    )
  

  plots[[i]] <- ggplot(visit_summary, aes(x = season, y = mean)) + 
    geom_col() + facet_grid(~type_visit) +
    ggtitle(paste("Average Visits by Season for", park_types[i]))
}

plots
## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

## 
## [[7]]

## 
## [[8]]

## 
## [[9]]

## 
## [[10]]

## 
## [[11]]

## 
## [[12]]

## 
## [[13]]

## 
## [[14]]

## 
## [[15]]

## 
## [[16]]

## 
## [[17]]

## 
## [[18]]